DNA Dynamics in A Water Drop 
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Due to its polyionic character the DNA double helix is stable and biologically active only in salty 
aqueous media where its charge is compensated by solvent counterions. Monovalent metal ions are 
ubiquitous in DNA environment and they are usually considered as the possible driving force of 
sequence- dependent modulations of DNA structure that make it recognizable by proteins. In an 
effort to directly examine this hypothesis, MD simulations of DNA in a water drop surrounded by 
vacuum were carried out, which relieves the requirement of charge neutrality. Surprisingly, with 
zero concentration of counterions a dodecamer DNA duplex appears metastable and its structure 
remains similar to that observed in experiments. 

PACS numbers: 87.14.Gg; 87.15. Aa; 87.15.He 
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It has been long recognized that, because of the 
polyionic nature of DNA, solvent counterions are re- 
quired for its stability, and that gross structural changes 
in DNA can be provoked by changing the concentration, 
the charge or the type of counterions |l| . More recently 
it has been proposed that free ions might also act as bio- 
logical regulators because, by binding to DNA, they can 
provoke conformational deformations recognized by spe- 
cific proteins [|[ ^, ^, |[ ^, |^. The most controversial 
is the role of the common monovalent cations Na"*" and 
K"*". They are ubiquitous in the DNA environment and 
can be readily available for any purpose. Until recently, 
they remained invisible in experimental DNA structures 
because it is difficult to detect them in water, and it 
has been suggested that they are perhaps responsible for 
the most widespread deformations of the double helix, 
namely, narrowing of the minor groove and bending. It 
is assumed that counterions are sequestered in the minor 
groove of some sequences, which brakes the symmetry of 
the repulsive electrostatic forces and provokes deforma- 
tions. This model is general and it easily explains other 
puzzling effects in DNA structure. 

The above hypothesis is supported by many recent 
studies. Penetration of monovalent cations into the mi- 
nor DNA groove has been confirmed by X-ray diffraction 
|, |, p^, N MR spectroscopy and MD sim- 

ulations y, However, it appears difhcult 
to find a discriminating set-up for testing the cause and 
consequence relationship between the ions and the DNA 
structure. All available evidences have more than one 
interpretation making this problem highly controversial 
p8| , [T9I , For example, correlations between ion posi- 
tions and the local groove width observed in MD cannot 
answer whether the ions perturb DNA or they just bind 
"opportunistically" in the sites of low potential near al- 
ready deformed double helix |Q. To clarify the issue 
of cause and effect one would have to remove solvent 
ions and check if the supposed counterion effects disap- 



pear with them. Unfortunately, a counterion-free DNA 
does not exist in nature whereas the most reliable compu- 
tational procedures presently employed require that the 
simulation cell that holds DNA carries zero net charge. 
Therefore, in both simulations and experiments the coun- 
terion effects cannot be completely eliminated. 

In an effort to directly address this question I have 
adapted the Particle Mesh Ewald algorithm ^ for 
modeling dynamics of DNA in a water drop surrounded 
by vacuum. Particle-mesh calculations with vacuum 
boundary conditions are long known in physics [Q, but, 
to my knowledge, they were never applied to chemical 
or biological systems. Free vacuum boundaries are in- 
tuitively most simple and they allow one to relieve the 
problems of charge neutrality and possible artifacts from 
interactions between periodical images. I describe here 
the first such "naive" all atom simulations of DNA in 
water, with unperturbed Coulomb electrostatics. It ap- 
pears that, with zero concentration of counterions, a do- 
decamer DNA duplex is metastable and its structure re- 
mains similar to that observed in experiments. 

The Dickerson-Drew dodecamer (CGCGAATTCGCG, 
in a canonical B-DNA conformation is sur- 
rounded by a spherical drop of 4000 TIP3P water 
molecules Initially, the drop had around 50 A in 

diameter and in dynamics it remained roughly spheri- 
cal. A rectangular unit cell is constructed around the 
drop with the minimal separation of 25 A between the 
water molecules and the cell sides. The cell is repli- 
cated periodically, which gives an infinite lattice of wa- 
ter drops with at least 50 A spacing between the clos- 
est neighbors. A shifted Coulomb's law is used, with 
U{rij) = ZiZj{l/rij - l/i?o) for < _Ro and [/(r^) = 
for longer distances, where Rq = 50 A, which eliminates 
any interactions between periodical images. Because this 
shifting does not affect the forces, the system behaves in 
dynamics as if surrounded by infinite vacuum. Within 
the drop the electrostatics are effectively evaluated with 
a cut-off of 50 A, which is larger than the DNA size and 
applies only to a small fraction water molecules at oppo- 
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TABLE I: Some structural parameters of standard and com- 
puted DNA conformations. Sequence avera ged helical param- 
eters were computed with program Curves |27| . All distances 
are in angstroms and angles in degrees. 





Xdisp 


Inclin 


Rise 


Twist 


RMSD-A" 


RMSD-B" 


A-DNA 


-5.4 


+ 19.1 


2.6 


32.7 


0.0 


6.2 


B-DNA 


-0.7 


-6.0 


3.4 


36.0 


6.2 


0.0 


Tjl" 


-2.2 


+6.9 


3.3 


33.8 


4.52 


2.36 


Tj2= 


-2.0 


+6.5 


3.2 


34.1 


4.57 


2.18 




-2.1 


+6.9 


3.3 


33.6 


4.46 


2.41 



"Heavy atom root mean square deviation from the corresponding 
canonical DNA form. 
''Water drop calculation without counterions. 
■^Water drop neutraUzed by 22 Na+ ions. 

'^Conventional PME calculation with counterions and periodical 
boundaries. 




FIG. 1: A snapshot from the last nanosecond of the trajec- 
tory with water oxygen positions shown by dots. 



site poles of the drop. Larger distances were also tested, 
but showed no noticeable difference. The Van-der-Waals 
interactions are computed with a conventional cut-off of 
9 A. 

The shifted Coulomb interactions are evaluated by us- 
ing the SPME method H]. The shifting is taken into 
account in the reciprocal sum by replacing the Fourier 
transform of 1/r with that of the shifted law. The re- 
sulting series is absolutely convergent regardless of the 
system charge. The values of energies and forces thus ob- 
tained have been verified by direct computation without 
cut-offs. Modifications in the direct sum are not neces- 
sary if we are only interested in the forces and not in the 



absolute energy values. It can be shown that this sim- 
ple procedure is essentially equivalent to mathematically 
more complex earlier formulations of the Ewald method 
for cluster calculations In dynamics, the shape 

of the drop fluctuates and the size of the unit cell is ad- 
justed accordingly at every time step. Evaporated water 
molecules are detected and excluded from calculations to 
prevent explosion in the cell size. After every 50 ps the 
calculation is stopped and water molecules that have left 
the drop are re-introduced with zero velocities by scatter- 
ing them randomly near the surface of the drop. The rate 
of evaporation was around 84 mol/ns, that is on average 
four molecules had to be re-introduced at each stop. 

All calculations were carried out with AMBER98 pa- 
rameters [|9[ |30| by using the ICMD method [|l[ || |3| 
with the time step of 0.01 ps. Two control simulations of 
the same dodecamer have been carried out for compari- 
son. In the first, 22 Na"*" ions were included in the same 
water drop. The second represents a conventional cal- 
culation in a rectangular unit cell including DNA, 3901 
water molecules, and 22 Na+ ions by the original SPME 
method with periodical boundary conditions. The tra- 
jectories were continued to 5 ns. 

Figure |^ shows a snapshot of the charged drop from 
the last nanosecond of the trajectory. The system ap- 
pears metastable in the nanosecond time range. The 
DNA molecule shows no signs of denaturation. The wa- 
ter media remains continuous without internal bubbles. 
The surface of the drop is covered by spiky "protuber- 
ances" which reduce with increased drop size and disap- 
pear if DNA is neutralized by counterions. The DNA 
structures averaged over the last 2.5 ns of the three tra- 
jectories are characterized in Table |. These data indicate 
that with zero counterion concentration this dodecamer 
DNA molecule remains in B-form. It is understood that 
with increased DNA length it should have exploded be- 
cause otherwise the electrostatic energy would go to in- 
finity. However, no such trend is seen in Table |, which 
indicates that critical lengths of catastrophic deforma- 
tions are much larger whereas the dodecamer B-form is 
metastable in water even if it is charged. 

The scales of differences between the DNA conforma- 
tions in the three trajectories are well characterized by 
the data in Table | and they were always far from statis- 
tically significant. The corresponding experimental data 
are available only for the twist, namely, its value drops by 
0.1° when the NaCl concentration is reduced from 0.3M 
to 0.05M |||. A shift by 0.1° in the average twist is too 
small to be detected in a 5 ns simulation because, for 
this relatively small molecule, its fluctuations between 
consecutive 1 ns averaged structures can reach 1.0°. The 
absolute twist value is lower than in experiment, which is 
a known general feature earlier discussed in the literature 

The most famous feature of this DNA molecule is the 
middle A ATT fragment. It is long known from experi- 
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FIG. 2: (a) The time evolution of the profile of the minor 
groove in Tjl. The surface plots are formed by time-averaged 
successive minor groove profiles, with that on the front face 
corresponding to the final DNA conformation. The groove 
width is evaluated by using space traces of C5' atoms 
Its value is given in angstroms and the corresponding canoni- 
cal B-DNA level of 7.7 A is marked by the thin straight lines 
on the faces of the box. Note that the groove width can be 
measured only starting from the third base pair from both 
termini. 

(b) The minor groove profile averaged over the last nanosec- 
ond of Tjl. The central trace represents the average groove 
width with rms fiuctuations shown as error bars. The upper 
and lower solid traces show the maximal and minimal values, 
respectively. The dotted trace exhibits the profile observed 
in the experimental X-ray structure [|^. The canonical B- 
DNA groove width is marked by the horizontal dashed line. 
Despite the narrowing, the groove width remains larger than 
the canonical value , which corresponds to the lower aver- 
age twist. 



ments that the minor DNA groove always narrows in this 
and some similar sequences, called A-tracts, and widens 
outside of them. Figure |^ exhibits dynamics of the minor 
groove in Tjl and its average profile over the last nanosec- 
ond. It has a characteristic waving shape with a narrow- 
ing in the middle. The amplitude of this modulation is 
similar to that in the experimental X-ray structure. The 
minimal width is 1.5 A larger than experimental value, 
which is probably linked mechanically to the lower av- 
erage twist. Similar results were obtained for the other 
two trajectories and they are close to earlier reported 
simulation studies carried out with non-zero counterion 
concentrations 



The sequence-dependent groove-width modulations in 
DNA are well established experimentally and, in the re- 
cent years, they have been proposed to result from in- 
teractions with bound monovalent metal ions commonly 
undetectable in X-ray crystal maps ||, |[ |l^ . The present 
results evidence that it is not the case, supporting recent 
conclusions of different groups |l9[ ^ . They explain 
also why groove modulations and intrinsic DNA bending 
could be reproduced in MD simulations with simplified 
treatment of electrostatic interactions that ignored spe- 
cific counterion effects ^ p8[ . 

According to the counterion condensation theory, DNA 
in aqueous environment should be always covered by a 
shell of counterions and its charge should be compen- 
sated by around 75% regardless of the bulk ion concen- 
trations |3^. The results presented here do not contra- 
dict this theory but they are somewhat at odds with an 
implicit assumption that the counterion cloud is critical 
for the native DNA structure. This was surprising, at 
least for the author, and suggests that we are still far 
from complete understanding of interactions that control 
the DNA structure. The correlations observed in ear- 
lier MD simulations pTf apparently are due to binding 
of counterions in sites of low potential near an already 
narrowed minor groove, therefore, these interactions are 
structure-specific rather than sequence-specific, and they 
cannot be the driving force of the corresponding DNA 
deformations. 
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